---
title: "Data Science Transport – Second Assignment – Group 12"
output:
flexdashboard::flex_dashboard:
orientation: columns
vertical_layout: fill
source_code: embed
---
```{r setup, include=FALSE}
library(flexdashboard)
library(gtfsrouter)
library(tidyverse)
library(tidytransit)
library(sf)
library(tmap)
library(units)
library(RColorBrewer)
tmap_mode("view")
```
Intro {data-icon="fa-hourglass-half"}
=====================================
Column {data-width=100}
-------------------------------------
### Last Subject
[Last time](https://bernd886.github.io/data-science-transport-mid-assign-wise-2021/) we took a look at the GTFS data of the VBB and the output data of car-free zone MATSim simulation.
We stuck to the data sets. This time you will see plots showing results of :
- **GTFS:** routing on the whole pt network (last time only BVG and S-Bahn)
- **MATSim:** emissions distribution
For the raw files see the [Github page](https://github.com/bernd886/data-science-transport-final-assign-wise-2021).
Column {data-width=300}
-------------------------------------
### Some of Last Time Plots
Isochrones {data-icon="fa-expand-arrows-alt"}
=====================================
```{r, include = FALSE}
##############################################################
#
# READ GTFS DATA
#
##############################################################
# set work directions
setwd_gtfs <- function(){setwd("~/Documents/Uni/Master/DataScienceTransport/data/vbb-gtfs")}
setwd_data <- function(){setwd("~/Documents/Uni/Master/DataScienceTransport/data")}
setwd_work <- function(){setwd("~/Documents/Uni/Master/DataScienceTransport/assignment_2")}
setwd_work
# read gtfs data for monday
file <- file.path("~/Documents/Uni/Master/DataScienceTransport/data/vbb-gtfs/2020-12_2020-12-28.zip")
gtfs <- extract_gtfs(file) %>% gtfs_timetable(day = 2)
##############################################################
#
# SET TIMES
#
##############################################################
start_time <- 7 * 3600 + 1200
end_time <- 8 * 3600
##############################################################
#
# CREATE STOPS SF OBJECT
#
##############################################################
stops <- st_as_sf(gtfs$stops,
coords = c("stop_lon", "stop_lat"),
crs = 4326) %>%
st_transform(25833)
##############################################################
#
# SHAPE DISTRICTS NEW (+ area)
#
##############################################################
setwd_data()
shape_districts_new <- read_sf(dsn = "LOR_SHP_2019-1", layer = "Planungsraum_EPSG_25833")
setwd_work()
shape_districts_new <- shape_districts_new %>%
group_by(BEZIRK) %>%
summarise() %>%
filter(!is.na(BEZIRK)) %>%
rename(NAME = BEZIRK) %>%
mutate(AREA = st_area(geometry)) %>%
select(NAME, AREA, everything()) %>%
mutate(AREA = (AREA / 1000000) * as_units("km2"))
# setting crs of polygons
st_crs(shape_districts_new$geometry) <- 25833
shape_berlin <- st_union(shape_districts_new)
##############################################################
#
# SHAPE CENTER AREAS
#
##############################################################
# "Zentrentragender Stadtraum mit höchster / hoher Urbanität"
# of Zentrumsbereichskernen
# see page 39: https://www.stadtentwicklung.berlin.de/planen/stadtentwicklungsplanung/download/zentren/2011-07-31_StEP_Zentren3.pdf
# or page 45 (less detailed): https://www.stadtentwicklung.berlin.de/planen/stadtentwicklungsplanung/download/zentren/StEP_Zentren_2030.pdf
# recreated with QGis
shape_center <- read_sf(dsn = "shape_center_areas", layer = "center_areas") %>%
mutate(name = c("east", "west")) %>%
select(name)
shape_center_east <- shape_center %>% filter(name == "east")
shape_center_west <- shape_center %>% filter(name == "west")
##############################################################
#
# SPECIFIC SHAPES AND STOPS
#
##############################################################
stops_in_berlin <- stops %>%
mutate(inside_berlin = st_within( geometry, shape_berlin )) %>%
mutate(inside_berlin = !is.na( as.numeric( inside_berlin ))) %>%
filter(inside_berlin == TRUE)
##############################################################
#
# CALCULATE ISOCHRONES
#
##############################################################
# # the following code calculates the isochrones (inefficent, ~ 15h)
# # instead of running the code, we read in the pre-calculated file
# stops_ic_area <- vector(mode = "double")
#
# # create isochrone areas for stops in 50 minutes
# for (stop_name in stops$stop_name){
#
# tryCatch( {
# ic_area <- gtfs_isochrone (gtfs,
# from = stop_name,
# #from_is_id = TRUE,
# start_time = start_time,
# end_time = end_time)$hull$area
# if(is.null(ic_area)) {
# stops_ic_area <<- rbind(stops_ic_area, 0)
# print(paste(stop_name, ": ", ic_area, "!!!!!!!!!!"))
# } else {
# stops_ic_area <<- rbind(stops_ic_area, ic_area)
# print(paste(stop_name, ": ", ic_area))
# }
# },
# error = function(e) {
# stops_ic_area <<- rbind(stops_ic_area, 0)
# print(paste("ERROR!!!", stop_name))
# }
# )
# }
#
# ##############################################################
# #
# # CLEANING
# #
# ##############################################################
#
# # merge and clean
# # https://r-spatial.github.io/sf/reference/bind.html
# # https://cran.r-project.org/web/packages/units/vignettes/units.html
# rownames(stops_ic_area) <- NULL
# stops_area <- st_sf(data.frame(stops, stops_ic_area / 1000000)) %>%
# rename(ic_area = stops_ic_area.1e.06,
# id = stop_id,
# name = stop_name,
# parent = parent_station) %>%
# select(id, name, parent, ic_area) %>%
# mutate(ic_area = ic_area * as_units("km2"))
#
# # save
# # https://r-spatial.github.io/sf/reference/st_write.html
# st_write(stops_area, "output_stops_ic_area.shp")
stops_area <- st_read("output_stops_ic_area.shp")
# more cleaning for plot
# https://dplyr.tidyverse.org/reference/distinct.html
stops_area =
stops_area %>%
select(name, ic_area) %>%
distinct(name, .keep_all = TRUE)
stops_area_berlin <- stops_area %>%
mutate(inside_berlin = st_within( geometry, shape_berlin )) %>%
mutate(inside_berlin = !is.na( as.numeric( inside_berlin ))) %>%
filter(inside_berlin == TRUE) %>%
select(-inside_berlin) %>%
mutate(id = paste(name, ": ", round(ic_area)))
```
Column {data-width=100}
-------------------------------------
### Approach {data-height=200}
The [gtfsrouter](https://atfutures.github.io/gtfs-router/) package allows us to calculate all stations reachable within a specified time period from a nominated station ([isochrones](https://atfutures.github.io/gtfs-router/reference/gtfs_isochrone.html)). We use the hull polygon as an indicator for the city area reachable.
See below for an isochrone plot of a 30min Monday work–home trip starting from Helmholtzstraße at 18:00.
We tried to simulate a home–work trip in the morning rush hour arriving at 08:00. Assuming a five minute walk at beginning and end of the trip, we limited the GTFS data to Monday (2021-01-18) and calculated the isochrones for a 40min time range. By that we wanted to get all reachable stations by local transport (< 60min, [PBefG §8 (1)](https://www.gesetze-im-internet.de/pbefg/__8.html)) and rate the station inside the network (reachability of possible work places).
**Some critique:** Since we've done the routing with all transport agencies and a large time range, the hull area is even more dominated by reachable fanning out long-distance transport. Additionally the area is not weighted by any form of utility.
### 30min IC Home Trip (18:00) from Helmholtzstr. {data-height=100}
```{r}
##############################################################
#
# Helmholtzstr. ISOCHRONE
#
##############################################################
ic_einstein <- gtfs_isochrone(gtfs,
from = "Berlin, Helmholtzstr.",
start_time = 18 * 3600,
end_time = 18 * 3600 + 1800)
tm_basemap(leaflet::providers$OpenStreetMap.DE) +
tm_shape(ic_einstein$hull) +
tm_polygons(col = "red",
alpha = 0.2,
border.col = "red") +
tm_shape(ic_einstein$routes) +
tm_lines() +
tm_shape(ic_einstein$end_points) +
tm_dots(col = "red") +
tm_shape(ic_einstein$start_point) +
tm_dots(col = "green")
```
Column {data-width=300}
-------------------------------------
### Area Size of the Hull Enclosing the Routed Points
```{r}
##############################################################
#
# PLOT
#
##############################################################
tm_shape(shape_districts_new) +
tm_polygons(alpha = 0,
popup.vars = c("area" = "AREA")) +
tm_shape(stops_area_berlin) +
tm_dots(col = "ic_area",
id = "name",
popup.vars = c("area" = "ic_area"),
size = 0.07,
border.lwd = 0.3,
legend.hist = TRUE,
n = 15,
title = "isochrone area [km^2]") +
tm_view(bbox = shape_center)
```
Travel Times to Center {data-icon="fa-stopwatch"}
=====================================
```{r, include = FALSE}
##############################################################
#
# READ GTFS DATA
#
##############################################################
# now we work with tidytransit
# calculation of shortest tt from all station to specific ones is more convinent
setwd_gtfs()
gtfs <- read_gtfs("2020-12_2020-12-28.zip")
setwd_work()
# http://tidytransit.r-transit.org/reference/filter_stop_times.html
stop_times_filtered <- filter_stop_times(gtfs, "2021-01-18", "06:00:00", "07:55:00")
##############################################################
#
# GET STOPS
#
##############################################################
stops <- st_as_sf(gtfs$stops, coords = c("stop_lon", "stop_lat"), crs = 4326) %>%
st_transform(25833) %>%
select(stop_name) %>%
rename(name = stop_name) %>%
distinct(name)
stops_berlin <- stops %>%
mutate(inside_berlin = st_within( geometry, shape_berlin )) %>%
mutate(inside_berlin = !is.na( as.numeric( inside_berlin ))) %>%
filter(inside_berlin == TRUE) %>%
select(name)
stops_center <- stops %>%
mutate(inside_center = st_within( geometry, shape_center )) %>%
mutate(inside_center = !is.na( as.numeric( inside_center ))) %>%
filter(inside_center == TRUE) %>%
select(name)
stops_center_east <- stops %>%
mutate(inside_center_east = st_within( geometry, shape_center_east )) %>%
mutate(inside_center_east = !is.na( as.numeric( inside_center_east ))) %>%
filter(inside_center_east == TRUE) %>%
select(name)
stops_center_west <- stops %>%
mutate(inside_center_west = st_within( geometry, shape_center_west )) %>%
mutate(inside_center_west = !is.na( as.numeric( inside_center_west ))) %>%
filter(inside_center_west == TRUE) %>%
select(name)
##############################################################
#
# TT calculation
#
##############################################################
# what are the tt to the center areas?
# according to Nahverkehrsplan Berlin 2019-2023: ANlage 1 - Monitoringbericht (p. 12)
# standard: tt_max = 3600, n_transfer_max = 2, n_realise_stations = 0.95
tt <- travel_times(
stop_times_filtered,
stops_center$name,
time_range = 5400,
arrival = TRUE,
max_transfers = 2,
# max_departure_time = NULL,
return_coords = TRUE,
return_DT = FALSE
)
# clean it for plot
tt <- tt %>%
rename(from = from_stop_name,
to = to_stop_name,
tt = travel_time,
departure = journey_departure_time,
arrival = journey_arrival_time
) %>%
select(-c(from_stop_id, to_stop_id, to_stop_lat, to_stop_lon)) %>%
st_as_sf(coords = c("from_stop_lon", "from_stop_lat"),
crs = 4326) %>%
st_transform(25833) %>%
mutate(tt = set_units(round(tt/60, 2), "min"))
```
Column {data-width=100}
-------------------------------------
### Approach {data-height=400}
The [local transport plan](https://www.berlin.de/sen/uvk/verkehr/verkehrsplanung/oeffentlicher-personennahverkehr/nahverkehrsplan/) (p. 106) sets targets for the connectivity standards. Different categories of center areas (see [StEP](https://stadtentwicklung.berlin.de/planen/stadtentwicklungsplanung/de/zentren/zentren2030/index.shtml), p. 45) have to be reachable within a certain time and with a maximum number of transfer. This should hold for 95% of the stations.
Based on the GTFS data, we tried to recreate the result of the monitoring ([NVP Anlage 1](https://www.berlin.de/sen/uvk/_assets/verkehr/verkehrsplanung/oeffentlicher-personennahverkehr/nahverkehrsplan/broschure_nvp_2019_anlage_1.pdf), p. 12), mentioning a degree of fulfillment of 99,7% for the central areas.
* destination:
+ City West (Zoo/ Kurfürstendamm)
+ Mitte (Potsdamer Platz/ Alexanderplatz)
* max. tt: 60min
* max. transfers: 2
From a more detailed [illustration of the area](https://www.stadtentwicklung.berlin.de/planen/stadtentwicklungsplanung/download/zentren/2011-07-31_StEP_Zentren3.pdf) (p. 39), we created a shape file enclosing the associated stations. The [tidytransit](https://tidytransit.r-transit.org/) package let us [calculate](https://tidytransit.r-transit.org/reference/travel_times.html) the shortest travel time for all stations to any of a specified set of stations. For that `arrival` had to be set to `TRUE`.
We wanted to simulate a home–work trip in the morning rush hour arriving at 08:00. Assuming five minute walking at beginning and end of the trip, we limited the GTFS data to Monday (2021-01-18) and calculated the shortest travel times in a 90min time range ending at 07:55.
### Percent of All Berlin Stations Fullfill the Connectivity Standard According to the GTFS Data {data-height=200}
```{r}
##############################################################
#
# DEGREE OF FULLFILMENT
#
##############################################################
n_of_stations <- tt %>%
mutate(inside_berlin = st_within( geometry, shape_berlin )) %>%
mutate(inside_berlin = !is.na( as.numeric( inside_berlin ))) %>%
filter(inside_berlin == TRUE) %>%
mutate(outside_center = st_within( geometry, shape_center )) %>%
mutate(outside_center = is.na( as.numeric( outside_center ))) %>%
filter(outside_center == TRUE) %>%
nrow()
n_of_stations_valid <- tt %>%
mutate(inside_berlin = st_within( geometry, shape_berlin )) %>%
mutate(inside_berlin = !is.na( as.numeric( inside_berlin ))) %>%
filter(inside_berlin == TRUE) %>%
mutate(outside_center = st_within( geometry, shape_center )) %>%
mutate(outside_center = is.na( as.numeric( outside_center ))) %>%
filter(outside_center == TRUE) %>%
filter(tt <= 60 * as_units("min")) %>%
filter(transfers <= 2) %>%
nrow()
percent_stations_valid <- n_of_stations_valid / n_of_stations * 100
percent_stations_valid <- round(percent_stations_valid, 2)
valueBox(paste(percent_stations_valid, "%"), icon = "fa-crosshairs")
```
Column {data-width=300}
-------------------------------------
### Shortest Travel Time to One of the Stations Inside City West or Mitte
```{r}
##############################################################
#
# PLOT
#
##############################################################
# https://campus.datacamp.com/courses/visualizing-geospatial-data-in-r/raster-data-and-color?ex=9
rdylgn <- rev(brewer.pal(7, "RdYlGn"))
# https://leaflet-extras.github.io/leaflet-providers/preview/
# https://tlorusso.github.io/geodata_workshop/tmap_package
# https://www.rdocumentation.org/packages/tmap/versions/3.0/topics/tm_basemap
# https://rdrr.io/cran/tmap/man/tm_view.html
# https://leafletjs.com/reference-1.3.4.html#map-methods-for-modifying-map-state
tm_basemap(leaflet::providers$CartoDB.DarkMatter) +
tm_shape(shape_districts_new) +
tm_polygons(alpha = 0,
lwd = 1.5,
border.col = "white",
popup.vars = c("area" = "AREA")
) +
tm_shape(shape_center) +
tm_polygons(alpha = 0.2,
col = "red",
border.col = "red"
) +
tm_shape(tt) +
tm_dots(col = "tt",
style = "fixed",
breaks = c(0, 10, 20, 30, 40, 50, 60, 120),
labels = c("0 – 10", "10 – 20", "20 – 30", "30 – 40", "40 – 50", "50 – 60", "> 60"),
id = "from",
palette = rdylgn,
title = "traveltime [min]",
popup.vars = c("to" = "to",
"traveltime" = "tt",
"departure at" = "departure",
"arrival at" = "arrival",
"number of transfers" = "transfers")
) +
tm_view(bbox = shape_center)
```
Emission Analysis for Berlin {data-icon="fa-burn"}
=====================================
```{r, include = FALSE}
berlin_bezirke <- st_read("~/Documents/Uni/Master/DataScienceTransport/assignment_2/Hao/shp-bezirke/bezirke_berlin.shp")
berlin_emissions <-read_delim("~/Documents/Uni/Master/DataScienceTransport/assignment_2/Hao/berlin-v5.5-1pct.emissionsgrid_Berlin_PlanA.csv",
delim="\t",
# sep = "\t",
locale=locale(decimal_mark = "."),
col_types = cols(
x = col_double(),
y = col_double()
))
map = read_sf("~/Documents/Uni/Master/DataScienceTransport/assignment_2/Hao/shp-bezirke/bezirke_berlin.shp")
berlin_emissions_sf <- st_as_sf(berlin_emissions, coords = c('x', 'y'), crs = st_crs(map))
berlin_emissions_bezirke <- berlin_emissions_sf %>% mutate(
intersection = as.integer(st_intersects(geometry, map))
, area = if_else(is.na(intersection), '', map$Name[intersection])
)
# ----
joined <- st_join(berlin_bezirke, berlin_emissions_bezirke)
```
Column {data-width=100}
-------------------------------------
### Emission Level per Districts for Pollutants eg. NH3 & PM2.5 {data-height=400}
```{r}
# Plot1: PM2.5+NH3 ----
joined_count_1 <- joined %>%
group_by(Name) %>%
summarise(sum=sum(PM2_5_non_exhaust)) %>%
#summarise(N2O_sum=sum(N2O)) %>%
#summarise(sum=sum(NH3)) %>%
ungroup() %>%
mutate(col_value=sum/n()) %>%
#mutate(col_value=1/sum*n())
rename("Fine Particulate Matter_average" = "col_value")
joined_count_2 <- joined %>%
group_by(Name) %>%
#summarise(sum=sum(PM2_5_non_exhaust)) %>%
#summarise(N2O_sum=sum(N2O)) %>%
summarise(NH3_sum=sum(NH3)) %>%
ungroup() %>%
mutate(col_value_2=NH3_sum/n()) %>%
#mutate(col_value=1/sum*n())
rename("NH3(Ammonia)_average" = "col_value_2")
tmap_mode("view")
tm_shape(joined_count_2) +
# tm_borders() +
tm_polygons(col="NH3(Ammonia)_average", palette = "Blues") + #Reds RdBu
tm_shape(joined_count_1) +
# tm_dots(size="size_value", col="col_value")
tm_bubbles(size="sum", col="Fine Particulate Matter_average")
```
### {data-height=100}
**Size of the Dots:** Sum of Fine Particulate Matter for this district
**Color of the Dots:** Average value of Fine Particulate Matter for this district
**Color of the Districts:** Average value of NH3 for this district
Column {data-width=100}
-------------------------------------
### The Local Distribution of the Pollutant N2O in Berlin {data-height=400}
```{r}
# Plot2 N2O ----
joined_count <- joined %>%
group_by(Name) %>%
#summarise(sum=sum(PM2_5_non_exhaust)) %>%
summarise(sum=sum(N2O)) %>%
#summarise(sum=sum(NH3)) %>%
ungroup() %>%
#mutate(col_value=sum/n())
mutate(col_value=1/sum*n())
berlin_emissions_bezirke_2 <- berlin_emissions_bezirke %>%
rename("N2O(Nitrous Oxide)" = "N2O")
tmap_mode("plot") #tmap_mode("view")
tm_shape(joined_count) +
tm_borders() +
# tm_polygons() +
# tm_polygons(col="col_value") +
# tm_dots(size="size_value", col="col_value")
tm_shape(berlin_emissions_bezirke_2) +
tm_dots(size=0.01, col="N2O(Nitrous Oxide)", border.lwd=NA)
# tm_dots(size=0.001, col="N2O", alpha=0.1)
```
### {data-height=100}
**Dark Red:** high value of this pollutant.
**Light Yellow:** low value of this pollutant.
Problem: Emissions of Charlottenburg
=> implement “Car-free Zone” transportation policy
Routing for Vehicles {data-icon="fa-car"}
=====================================
Column {data-width=100}
-------------------------------------
### Background
Car-free Area Policy: Only the residents can use their cars inside of the Car-free Area. Non-residents like workers are not allowed to drive inside of the Car-free Area.
Scenario: Please consider a worker who works inside our study area (Car-free Area). The Red Dot stands for his/her/their office location. This Person wants to go home after work.
In the Base case: He/She can take any transportation mode for going back home. For example, driving by his/her car to get back home.
In the Policy Case (Car-free Area Case): He/She can not go back home by driving anymore. The alternatives can only be public transportation in most cases.
Column {data-width=300}
-------------------------------------
```{r, include = FALSE}
library(sfnetworks)
library(tidygraph)
library(TSP)
tmap_mode("view") # tmap_mode("plot")
berlin_network <- sf::st_read(
dsn = "https://download.geofabrik.de/europe/germany/berlin-latest.osm.pbf",
# dsn = "/Users/haowu/Workspace/R/DataScience_HA2/Berlin.osm.pbf",
# dsn = "https://download.geofabrik.de/europe/great-britain/england/greater-london-latest.osm.pbf",
layer = "lines",
# query = "SELECT * FROM lines WHERE (highway != 'NA') AND (highway != 'unclassified')",
query = "SELECT * FROM lines WHERE (highway = 'primary') OR (highway = 'secondary') AND (highway != 'NA') AND (highway != 'unclassified')",
# query = "SELECT * FROM lines WHERE (highway = 'primary') OR (highway = 'secondary') OR (highway = 'residential') AND (highway != 'NA') AND (highway != 'unclassified')",
stringsAsFactors = FALSE
)
```
### Routing for Vehicles
```{r}
# net = as_sfnetwork(h2, directed = FALSE) %>%
net = as_sfnetwork(berlin_network, directed = FALSE) %>%
st_transform(25833) %>%
activate("edges") %>%
mutate(weight = edge_length())
# ----
# How many edge types are there?
types = net %>%
activate("edges") %>%
pull(highway) %>%
# filter(highway!=NA) %>%
# filter(!is.na(highway)) %>%
unique()
# Randomly define a driving speed in m/s for each type.
# With values between 18 and 30 km/hr.
set.seed(1)
speeds = runif(length(types), 18 * 1000 / 60 / 60, 30 * 1000 / 60 / 60)
# Assign a speed to each edge based on its type.
# Calculate travel time for each edge based on that.
net = net %>%
activate("edges") %>%
group_by(highway) %>%
# mutate(speed = units::set_units(speeds[cur_group_id()], "m/s")) %>%
mutate(speed = units::set_units(speeds[cur_group_id()], "m/s")) %>%
mutate(time = weight / speed) %>%
ungroup()
# ----
net = activate(net, "nodes")
p = net %>%
st_geometry() %>%
st_combine() %>%
st_centroid()
iso = net %>%
# filter(node_distance_from(st_nearest_feature(p, net), weights = time) <= 600)
filter(node_distance_from(st_nearest_feature(p, net), weights = time) <= 2400)
iso_poly = iso %>%
st_geometry() %>%
st_combine() %>%
st_convex_hull()
###--
tm_shape(st_as_sf(iso, "edges")) + tm_lines(col="highway", lwd = 5) +
tm_shape(st_as_sf(iso, "nodes")) + tm_dots(col="grey") +
tm_shape(iso_poly) + tm_polygons(col="black", alpha=0.1) +
tm_shape(p) + tm_dots(col = "red", size=0.02)
```
Sources {data-icon="fa-external-link-alt"}
=====================================
### Sources and Usefull Links
* VBB GTFS data:
+ [Website](https://www.vbb.de/unsere-themen/vbbdigital/api-entwicklerinfos/datensaetze)
+ [used data set](http://transitfeeds.com/p/verkehrsverbund-berlin-brandenburg/213/20201228) on transitfeeds
* Berlin shape files: [Senatsverwaltung für Stadtentwicklung und Wohnen](https://web.archive.org/web/20190624123508/https://www.stadtentwicklung.berlin.de/planen/basisdaten_stadtentwicklung/lor/de/download.shtml)
* GTFS routing:
+ [tidytransit](https://tidytransit.r-transit.org/) package
+ [gtfs-router](https://atfutures.github.io/gtfs-router/) package
* helpful references:
+ sf reference: [bind](https://r-spatial.github.io/sf/reference/bind.html), [st_write](https://r-spatial.github.io/sf/reference/st_write.html), [st_read](https://r-spatial.github.io/sf/reference/st_read.html)
+ [Units of Measurement for R Vectors](https://cran.r-project.org/web/packages/units/vignettes/units.html)
+ tidyverse dplyr: [distinct](https://dplyr.tidyverse.org/reference/distinct.html)
* connectivity standards:
+ [Stadtentwicklungsplan Zentren 3](https://www.stadtentwicklung.berlin.de/planen/stadtentwicklungsplanung/download/zentren/2011-07-31_StEP_Zentren3.pdf) (2011)
+ [Stadtentwicklungsplan Zentren 2030](https://www.stadtentwicklung.berlin.de/planen/stadtentwicklungsplanung/download/zentren/StEP_Zentren_2030.pdf) (2019)
+ [Nahverkehrsplan 2019-2023](https://www.berlin.de/sen/uvk/verkehr/verkehrsplanung/oeffentlicher-personennahverkehr/nahverkehrsplan/) (2019)
* plotting:
+ datacamp: [Custom palette in tmap](https://campus.datacamp.com/courses/visualizing-geospatial-data-in-r/raster-data-and-color?ex=9)
+ [Leaflet-providers preview](https://leaflet-extras.github.io/leaflet-providers/preview/)
+ tlorusso: [tmap](https://tlorusso.github.io/geodata_workshop/tmap_package)
+ rdocumentation: [tm_basemap](https://www.rdocumentation.org/packages/tmap/versions/3.0/topics/tm_basemap)
+ tm_view: [Options for the interactive tmap viewer](https://rdrr.io/cran/tmap/man/tm_view.html)
